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We present a model of anomalous diffusion consisting of an ensemble of particles undergoing 
homogeneous Brownian motion except for confinement by randomly placed reflecting boundaries. 
For power-law distributed compartment sizes, we calculate exact and asymptotic values of the 
ensemble averaged mean squared displacement and find that it increases subdiffusively, as either a 
power or the logarithm of time. Numerical simulations show that the probability density function 
of the displacement is non-Gaussian. We discuss the relevance of the model for the analysis of 
single-particle tracking experiments and its relation to other sources of subdiffusion. In particular 
we discuss an intimate connection with diffusion on percolation processes. 
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In the study of diffusion in physical systems it is well 
known that a variety of processes lead to diffusion that 
deviates from pure Brownian motion [1-6]. The most 
commonly studied signature of anomalous diffusion is the 
mean squared displacement (MSD), which often takes the 
form of a power law, 

( 1 ) 

where x{t) is the displacement of a particle at time t. 
For Brownian motion /? = 1, and for subdiffusive motion 
P < 1. Subdiffusive systems may be successfully mod¬ 
eled by various random processes such as continuous time 
random walk (CTRW) [7-9], fractional Brownian mo¬ 
tion [10, 11] (FBM), diffusion on percolation [12] (DOP), 
Lorentz models [5], or by a combination of these [13, 14]. 
The problem of distinguishing which processes contribute 
to or modify anomalous diffusion is an active area of re¬ 
search [15-18]. In particular, one may need to model 
complicating factors such as confinement. Confinement 
plays an important role in modifying or attenuating 
the subdiffusion manifested by the power-law behavior 
in (1) [19-28]. But, in this paper we take a different view 
of the relation of confinement to subdiffusion: We show 
that random confinement may in fact be the sole cause 
of observed subdiffusion. While the effect is quite gen¬ 
eral, the model and analysis are motivated specifically by 
single particle tracking (SPT) experiments, in particular 
studies of the biophysics of cells [29-32]. We present a 
scenario in which disordered confinement produces quan¬ 
titative signatures of subdiffusion in a typical analysis of 
SPT data. Given that heterogeneous confining bound¬ 
aries are often observed in biophysics, it follows that in 
experiments searching for contributions to the subdiffu¬ 
sion exponent /3, confinement may not be discarded a 
priori as a candidate. 

Experimental scenario — Consider a number of Brown¬ 
ian particles uniformly distributed in a space that is par¬ 
titioned into compartments by a random arrangement 
of reflecting barriers. Apart from the presence of the 


boundaries, the motion is diffusive, with parameters ho¬ 
mogeneous is space and time. We collect an SPT tra¬ 
jectory (time series) for each of several particles sam¬ 
pled uniformly from the space. The trajectory consists of 
the displacement of each particle from its starting point 
recorded at a series of times. Typically, the trajecto¬ 
ries are first analyzed via the ensemble averaged mean 
squared displacement (EMSD), or the time-ensemble av¬ 
eraged mean squared displacement (TEMSD), where the 
TEMSD consists, operationally, of first computing a time 
(sliding) averaged MSD (TMSD) for each trajectory and 
then averaging the result over the trajectories. For the 
minimal model introduced below, we find 

• The EMSD is unbounded and subdiffusive in the 
sense of (1) if the distribution of the linear size of 
the compartments has a sufficiently heavy-tail. 

• There is no weak ergodicity breaking (if the sys¬ 
tem is in equilibrium): the EMSD is equal to the 
TEMSD. 

• The TMSD (ie for a single trajectory) tends to a 
constant. That is, the observables of a trajectory 
are not self-averaging over disorder. 

The model — We choose a minimal model that cap¬ 
tures the essential features and displays asymptotically 
subdiffusive motion. We refer to the particular model 
presented below as the random scale-free confinement 
model (RSFC). We state the results for the MSD be¬ 
fore giving the detailed calculations. In this paper, we 
treat in detail only one spatial dimension, stating some 
results for higher dimension at the end, but leaving de¬ 
tails to a subsequent paper. The probability density for 
the displacement of a particle diffusing on a line segment 
of length r with reflecting boundary conditions obeys 

dtu{x;t) = dlu{x;t), (2) 

with dxu{r/2;t) = dxu{—r/2-,t) — 0. We assume that 
the probability for a particle to be found in a segment of 
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Figure 1. EMSD for random confinement model RSFC with 
ro = 1 showing bounded motion (a = 4, lower curve), power- 
law (a = 1, upper curve), and logarithmic (a = 2, middle 
curve) growth. Each of the three curves has points from exact 
solutions ((12) and (15)) and simulations. Green curve is 
asymptotic solution (14) for a = 1. Simulations are a rescaled 
random walk with ro = 100 in lattice units. Averages are over 
1.8 X 10® trajectories, and 10® steps per trajectory, using the 
Mersenne Twister random number generator. 



length r is Pareto-distributed 


P{r) = Pa{r) = 


0 for r < To 


arn r 


Oirr^ — a—l 


for r > To 


, for a > 0. 


( 3 ) 


We find that the asymptotic EMSD is, for a 2 

(^)=K^+K 2 t"-^, (4) 


so that the EMSD either decays to a constant for a > 2 
or grows without bound for a < 2. For a = 2 we find 


{x^{t)) = Ks + K 4 ln{t). (5) 


In (4) and (5) angle brackets denote averaging over parti¬ 
cle trajectories, the tilde denotes averaging over the dis¬ 
order (ie over r), and Ki, K 2 , K^, K 4 are constants that 
depend on dimension, and the initial distribution of par¬ 
ticles. Examples of curves for three values of a showing 
both bounded and unbounded EMSD are shown in Fig. 1. 


Comparison with other sources of subdijfusion — Our 
random confinement model does not manifest weak er- 
godicity breaking [33] , a characteristic whose presence is 
known to imply dramatic differences between EMSD and 
TEMSD observed in SPT experiments [34-39]. On the 
contrary, the model is in equilibrium at t = 0, and the 
increments of x{t) are stationary so that the TEMSD and 
EMSD are equal at all times. Thus, the subdiffusion is 
due to correlated increments, as is the case for both DOP 
and FBM. Our simulations show that the PDF of the 


displacement for large times is non-Gaussian as shown in 
Fig. 2, a feature shared by DOP, and CTRW, but not 
FBM. RSFC and DOP (but not FBM and CTRW) are 
examples of diffusion on disordered, or more broadly, het¬ 
erogeneous media. A model of this type that is closely 
related to RSFC considers domains of random scale-free 
size and diffusivity [40] . Despite this similarity, the latter 
model has deeper similarities to other models of heteroge¬ 
neous transport coefficients [6, 41-47]. When these mod¬ 
els do show subdiffusion, it arises from non-stationary 
increments, and they are thus more closely related to 
CTRW. 

RSFC and percolation — It is interesting to compare 
RSFC more closely to DOP. We distinguish two cases: 
The first is diffusion of a particle starting at a randomly 
chosen site on a percolation process at the critical thresh¬ 
old [12] (DOP I). In the second case the initial site is 
a random site on the critical infinite cluster (DOP II). 
Both DOP I and DOP II show subdiffusive EMSD, but 
with differing values of the anomalous diffusion exponent 
P in (1) [12, 48, 49]. For DOP II, subdiffusion is due 
to the fractal properties of the infinite cluster. How¬ 
ever, because the volume fraction of the critical infinite 
cluster is zero, it gives no weight in DOP I where the 
initial position may be any point on the lattice. As is 
the case with the pure random confinement of RSFC, 
in DOP I every particle is confined to a region of finite 
size. Thus, for both RSFC and DOP I, the time av¬ 
eraged mean squared displacement (TMSD) —the time 
average of a single trajectory— tends to a constant at 
long times. In DOP I, the size of the confinement regions 
(the finite clusters) has a heavy-tailed, power-law distri¬ 
bution. Thus, as for RSFC, averaging over an ensemble 
of uniformly distributed particles gives an unbounded, 
subdiffusive MSD, with an additional contribution from 
the fractal structure of large clusters. In other words. 
The fact that fi differs between (I) the walk on all clus¬ 
ters and (II) the walk on the infinite cluster, is due to 
random, scale-free confinement and RSFC represents an 
abstraction of this phenomenon. 

We address a potential issue in simulations and exper¬ 
iment. In CTRW, a finite number of trajectories show a 
subdiffusive EMSD no matter how long their duration, 
though this implies an ever increasing spatial domain. In 
DOP I and RSFC, the average over a finite number of 
trajectories will tend to a constant EMSD at long times, 
even for a < 2. But, simulations and SPT trajectories 
are typically limited to times shorter than the time T re¬ 
quired to explore the entire experimental domain. Since 
a heavy-tail (a < 2) implies that there are compartments 
whose area is of the order of the experimental domain, 
the time required to cross these domains will be of the 
same order as T, and thus power-law subdiffusion will 
be observed as shown in Fig. 1. Note that: In CTRW, 
no particle is trapped for an infinite time, yet the step 
rate decreases toward zero. In RSFC, the motion of every 
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Figure 2. Monte Carlo computation of PDF of displacement 
averaged over segment length r for parameters a = 1, ro = 
1, t = 100. Two lower curves are linear scale. Two upper 
curves are double log scale. Peaked curves; black linear scale 
and green log scale are equilibrium density (9) averaged over 
disorder. Flat-top curves: Violet linear scale and blue log 
scale, are average over disorder of (8) with xq — 0, ie non¬ 
equilibrium u{x,0) = S{x). Note that averaging over xq to 
get the equilibrium distribution increases the width of central 
region by a factor of two. 


particle is bounded for an infinite time, yet the EMSD 
grows without bound. 

Calculation of EMSD — We denote by [xo+y, t) the 
solution to (2) with r = 1 and initial condition u{x, 0) = 
6{x — xo)- We translate the origin to xq, and average 
over Xq to get the density for displacement from initial 
position for a trajectory sampled from the equilibrium 
distribution, which we denote by u{y,t). For y > 0 this 
average is 


u>{y,t) 


CD-y 

/ Uxg(xo + y,t)dxo 

J-l/2 


for 0 < ?/ < 1. 


may be written as an eigenfunction expansion 

OO 

u{x;t) = ^ ‘ cos(27rux) 

n—0 

sin(2^[n + 1/2])], 

with coefficients determined by the initial condition and 
orthogonality of the eigenfunctions. Using the initial con¬ 
dition u{x,0) = 6{x — Xo), we hnd 

OO 

m(x; Xq, t) = 1 + 2 ^ " * cos(27rnx) cos(27rnxo) 

n—1 

OO 

-1-2 (n.-i-i/2) ‘ sin(27r[n-b l/2]x) sin(27r[n-I-l/2]xo) 

n—0 

( 8 ) 


Evaluating the integrals obtained by substituting (8) into 
(6) is straightforward and results in 


u{y,t) =(1 - |y|) 


n—\ 


(1 - |y|)cos(7rn|2/|) 


sin(n7r|y|) 

Trn 


(9) 

The asymptotic shape of u(y, t) is reflected in the central 
portion of the disorder-averaged density as seen in Fig. 2. 
Inserting (9) into (7) we find 


1 1 

(x2(t;l)) = --— 

n=0 


>+l/2)^t 


(u +1/2)4. 


Setting t = 0 shows that 


1 ^ 1 1 

7r4 in -b 1/2)4 0 

n—0 ^ ' 


( 10 ) 


We now consider a segment of length r, rather than 1, 
and use (10) to write 


Because u(?/, t) is symmetric in y, we have 

r-l/S-EI 


r.2 °° 


u{y,t) = Uj;g{\y\+xo,t)dxo for - 1 < y < 1, 

J-l/2 

( 6 ) 

We denote by (x^(t;r)) the EMSD averaged over xq for 
a segment of length r. For the segment of unit length, 
this is given by 

(x^(t; 1)) = 2 / u{y,t)y'^dy. (7) 

Jo 

The EMSD averaged over the disorder is given by 


(x2(t)) = / (x^(t; r))P(r) dr. 

Jo 

Separating time and space variables, the solution to (2) 




n—0 


(n-b 1/2)4’ 


where we have absorbed the initial term into the sum in 
order to cancel two diverging quantities below. Averag¬ 
ing over the disorder we have 


1 


TW) =-i E 


1 


TT^ ^ (n + 1/2)^ 
n=0 V ' / ^ 

poo 

j (l-e-4-'(”+i/2)'*/’-')r2p(r)dr 


. OO 1 

TT 


TT ^fn-bl/2)‘ 

n=0 


/•OO _ 

/ (l — e“^) y~^^^P f27r(n -b tly \ dy. 

( 11 ) 
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Inserting the Pareto distribution (3) into (11), we have 


(xHt)) = 


27r4 ^ (n + 1/2)4' 


n—0 

^2 


_ ^'0 ^ 
“27r4 2^ (n + 1/2)4 

n—0 ^ ' ' 


Ein 


/ (1-e My 

Jo 


: — 2 47r^(n + l/2)^t 


( 12 ) 

where z = 47r^(n + and we have defined a gen¬ 

eralized entire exponential integral by 

Ein(s,z)=z“® f (l — e“^) y'’“^ dy. (13) 

Jo 


The series obtained by expanding the exponential[50] is 
easily seen to be analytic for s,z G C except for simple 
poles at s = —1, —2,.... In particular, (13) is analytic at 
s = 0, where Ein(0,z) = Ein(z), the usual entire expo¬ 
nential integral. This corresponds to the critical value of 
the power-law exponent a = 2. It can be shown that for 
a ^ 2, Ein(s, z) = 1/s — z~^j{s, z), where 7(5, z) is the 
lower incomplete gamma function. Large t corresponds 
to large z, so that, using lim^_,.oo 7(s,z) = r(s), we find 
the asymptotic form 


{xHt)) 


a-2 6 


- a7r-“-224-“r;) 



C{2 + a,l/2)t 



(14) 


where is the Hurwitz Zeta function. Eq. (14) 

shows that for a > 2, {x'^{t)) converges to a constant 

at long times, while for 1 < a < 2, {x'^{t)) shows sub¬ 
diffusion. For a = 2, we use Ein(0, z) = Ein(z) = 
ln(z) -I- 7 -I- r( 0 , z) and lim^_>oo r(s, z) = 0 to find the 
asymptotic solution 


(M(t)) 


7 -I- 2 In 



12a,C(4,l/2) 

■7r4 


-t- ln(t) 


(15) 

where 9sC(s, l/2)|s=4 « 10.9697. Equations (14) and 
(15) are the goals of the calculations and are shown in 
Fig. 1. 

Non-equilibrium initial conditions — Although we as¬ 
sumed equilibrium above, subdiffusion is also observed 
for other initial distributions. We performed calcula¬ 
tions analogous to those above for the initial distribution 
u{x, 0) = S{x), that is, each particle begins at the center 
of the confinement domain, and obtained the asymptotic 
form for a 7^ 2 


{x^t)} 


a-212 


- a7r-“2i-“r“r 



r]{a)t- 


( 16 ) 


where r]{x) is the Dirichlet eta function. For a = 2 we 
found 


{x^t)) 


rp f ln(t) 

6(2 


6C'(2) 

'TT^ 


-I- ln(7r) 


7 

2 



where ((x) is the Riemann zeta function. In two dimen¬ 
sions, ie diffusion on disks of random radius, we obtained 
the asymptotic forms, for a 7^ 2 


(pHt)} 


« 4-g 

a- 2 2 




where is the nth positive zero of the Bessel function 
of the first kind Ji^, and CmiP) is 


OAU) = E 


T t ■ \ ’ 

= 1 0 i^+i\ji^,n)jv,n 


and p{t) is the radial density, ie averaged over azimuthal 
angle. For a = 2, we found 

{7 - In(rg) - 165;3 Cji( 2) + ln(t)} , 

where dpCn{‘2) « 0.147342. 

Conclusion — We have demonstrated that heavy-tailed 
random scale-free confinement gives rise to a subdiffusive 
EMSD. We have discussed its relation to other sources of 
anomalous diffusion, in particular its presence in certain 
diffusive processes on percolation. A number of questions 
remain, including: Are there random potential fields that 
lead to RSFC ? What is nature of the propagator (PDF) 
of RSFC ? ... or of autocorrelation of the displacement? 
What is the anomalous exponent when RSFC is com¬ 
bined with other sources of anomalous diffusion ? 
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